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1 Introduction 



In order to study questions of chiral symmetry breaking [1-6] and universality 
[7,8], and also in the context of Liischer's fermion algorithm [9-11], one is 
interested in the eigenvalues of the gauge covariant Dirac operator (or 75 times 
the Dirac operator) which are close to the origin. One candidate method which 
comes to mind to determine eigenvalues numerically is a Lanczos algorithm 
[12]. Variants of this method have been used in lattice field theory for a long 
time, see e.g. Refs. [13,2] for staggered fermions and Refs. [14,3] for Wilson 
fermions. 

Presently there is renewed interest in the Lanczos method, for instance in 
Refs. [11,6], but the topics of these papers are such that the algorithmic fea- 
tures and difficulties connected with a Lanczos procedure were not addressed. 
The present paper complements this point. We perform a further algorithmic 
investigation of the Lanczos method in the context of lattice gauge theory 
which extends previous studies. We are also interested in a comparison of 
the partially converged Lanczos method (i.e. the case that just eigenvalues 
close to the origin are to be computed) with an accelerated conjugate gradient 
algorithm for low-lying eigenvalues [15,16]. 

The focus of the present paper will be on CuUum's and Willoughby's Lanczos 
method [17] applied to euclidean Wilson fermions in quenched and unquenched 
SU(2) gauge fields on lattices of volume ranging from 4"* to 16"*. The operators 
under consideration are 

g = 75(i) + m)/(8 + m) (1) 

and its square. [D -\- m) is the massive Dirac operator for Wilson fermions, 

[{D + m)il)]{x) = 
1 1 ^ 

At=i 

+ {l + -f^)U{x,x - ij)ij{x - fi)} . (2) 

Here k, = (2m -|- 8)~^ denotes the hopping parameter, x ± fi is the nearest 
neighbour site of x in it //-direction, and U{x^ x ± fi) is the gauge field on the 
link (x, X ± fi). In the present paper we impose periodic boundary conditions. 
The operator Q is hermitean, and it is normalized such that its eigenvalues 
are between -1 and 1. 

In some situations one does not have to know the sign of the eigenvalues of 
Q but only their absolute magnitude. This means that one can equally well 
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determine (a part of) the spectrum of . For example, one can probe the 
chiral hmit by means of the "pion norm" which depends on the eigenvalues 
of as discussed in Ref. [1]; or in Liischer's fermion algorithm [9,10] for an 
even number of flavours the small eigenvalues of can be used to correct 
for possible systematic errors in case that the polynomial approximation to 
the function 1/s is too poor at the lower end of the spectrum. If one is in 
such a situation, one can run the Lanczos procedure either with Q or with . 
We show that it is advantageous to diagonalize the unsquared operator which 
means that one could speed up computations as in Ref. [11]. 

Further points which will be discussed are given in the following overview. We 
concentrate on algorithmic aspects, physical questions will be addressed in a 
subsequent paper. In Sec. 2 we start by recalling CuUum's and Willoughby's 
Lanczos method. This method is reliable even on larger lattices, but its cost 
for the computation of a given fraction of the spectrum grows (approximately) 
with the square of the lattice volume. Issues of the convergence behaviour of 
the partially converged Lanczos method on lattices up to IG"* are described in 
Sec. 3. They reveal again the well-known fact that there is no black-box Lanc- 
zos routine, but one needs experience with the operator under consideration. 
In Sec. 4 we turn to the computation of complete spectra. It is shown that 
the convergence behaviour is closely related with the local spectral density. 
Complete spectra are determined on lattices up to 8"^ • 12 (and almost on a 
12"* lattice where 165884 of the 165888 eigenvalues were found). For configu- 
rations where all eigenvalues are computed, we give results for the fermionic 
determinants in Sec. 5. In this section we also quote values of determinants of 
quenched and unquenched staggered fermions whose spectra were determined 
in a previous publication [18]. Finally, we present spectral densities in Sec. 6, 
and we end with some conclusions and a comparison with the algorithm of 
Ref. [16]. 



2 Cullum's and Willoughby's Lanczos method 

The Lanczos procedure is a technique that can be used to solve large, sparse, 
symmetric or hermitean eigenproblems ^ [12]. The idea is to transform a given 
hermitean n X n matrix A into a similar symmetric tridiagonal matrix T = 
V~^AV with unitary and then T is diagonalized. The transformation of A 
can be performed iteratively. If one writes V = (^1,^2, • • • ,^^ra) with column 



^ In Refs. [13,14,3] the authors used a non-hermitean Lanczos method and come to 
the conclusion that it works on small lattices. 
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vectors Vi ("Lanczos vectors") and 



/ 



T 



/3i a2 



l-^n — l 



(3) 



then AV = VT is equivalent to 



Avi = aivi + /3iV2 , 

Av, = I3,_iv,_i + a,v, + I3,v,+i for z = 2, . . . , n - 1, (4) 

AVn = /3n-lVn-l + a^Vn ■ 

Given an initial — generally random — vector Vi and using the orthonormality 
among the Vi with respect to the canonical scalar product (•, •), one can 
determine iteratively from these equations: = {vi^Avi)^ jSf = [{vi^A^Vi) — 
— (with /3o = 0), and Vi^i = l3~^{Avi — — aiVi) as long as 

l3i ^ (otherwise the iteration is stopped). Note that T will be a real matrix 
also in case that A is complex. The sign of jSi is arbitrary. 

In exact arithmetic, the Lanczos iteration should finish after at most n steps 
and the last equation in (4) would be automatically fulfilled. For this case there 
exists a convergence theory for which we refer to the literature [12]. In practice, 
however, there are severe problems with a straightforward implementation of 
the Lanczos procedure [12,17]. These problems are caused by rounding errors 
and loss of orthogonality among the Lanczos vectors. ^ In principle the latter 
problem can be circumvented by storing all the Lanczos vectors and enforcing 
orthogonality among them by hand. However, then one is restricted to small 
lattices because of computer memory or I/O limitations. 

In CuUum's and Willoughby's proposal [17] one performs no reorthogonaliza- 
tion, and one continues the iteration (4) for an a priori unspecified count. In 
this way a sequence of j X j tridiagonal matrices T^^\ j = 1, 2, ... is generated. 
The diagonal elements of T^-'^ are a^, z = 1, . . . , j, and the off-diagonal entries 
are z = 1, . . . , j — 1. As a technical point we note that the a^, jSi and Vi^i 
will not be computed according to the formulas given above, but rather there 
exists a particular form of the recursion which has proven most stable, see [17] 
and [12, Algorithm 9.2.1 and remark on p. 492]. 



This loss of orthogonality is not necessarily due to the accumulation of round-off 
errors [12]. 
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The common belief is that generally the extremal eigenvalues of T^-'^ with 
increasing j are progressively better estimates of the extremal eigenvalues 
of A. Eventually, all the different eigenvalues of A will appear as eigenvalues 
of some r^-'^ with j generally larger than n. However, because of rounding 
errors and loss of orthogonality among the Lanczos vectors, there will also 
appear so-called "spurious" eigenvalues. These are eigenvalues of T^-'^ but not 
(approximate) eigenvalues of A. 

For the solution to the problem of identifying spurious eigenvalues and coping 
with their presence, CuUum and Willoughby give the following recipe [f7]. 
One compares the eigenvalues of T^-'^ with the eigenvalues of a matrix rj"*^ 
which equals T^-'^ with the first row and first column deleted. If a simple 
eigenvalue of T^-'^ is also an eigenvalue of T^''^ then this eigenvalue is spurious. 
(We remark that the authors of Refs. [13,14,2] work with T^-'"^) instead of 
rj''^ They identify approximate eigenvalues of A by having the property to be 
eigenvalues of both T^-'"^) and T^-'^ whereas spurious eigenvalues are different.) 

A problem which remains in the CuUum- Willoughby algorithm is that the 
multiplicities of the eigenvalues of T^-'^ do not refiect the correct multiplicities 
of the eigenvalues of A. However, this is neither a problem for the present 
study as we shall see, nor was it a problem in a previous study for staggered 
fermions [18]. 



3 Convergence behaviour of the Lanczos method for eigenvalues 
close to and furthest from the origin 

Henceforth the term Lanczos method refers to CuUum's and Willoughby's 
variant. The tridiagonal matrices T^^^ and r2"'^ were diagonalized by means of 
the Numerical Recipes routine "tqli" [19] which implements the QL algorithm 
with implicit shifts. Two eigenvalues were counted as different when they 
differed by more than I0~^°. This number is arbitrary but it is chosen such that 
it is small compared with the gaps in the spectra, and large compared with 
round-off errors. The computer program was checked for gauge covariance, 
and it was also verified that free spectra are obtained correctly, except for 
multiplicities. 

We turn first to the convergence behaviour of the Lanczos method for the 
operator Q defined in (I), and for eigenvalues A close to and furthest from the 
origin. The eigenvalues of T^-'^ were monitored as a function of j. Some results 
in individual gauge field configurations'* are collected in Table I. The entry 

The numerical data presented in this section were obtained from points in the 
/3-K-plane which are not in the immediate neighbourhood of a critical point. 
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Table 1 

Number of different eigenvalues (EVs) found after j iterations of the Lanczos pro- 
cedure for the operator Q. The examples are taken at k = 0.15, and (3 denotes the 
coupling constant of the SU(2) gauge field part of the action. 





124, ^ ^ 


2.4, quenched 


164, /3 = 


2.12, unquenched 


164, /3 = 


2.3, unquenched 


i 


EVs r(j') 


"good" EVs 


EVs r(j') 


"good" EVs 


EVs r(j') 


"good" EVs 


500 


500 


500 


500 


500 


500 


500 


1000 


992 


988 


1000 


998 


1000 


998 


2000 


1945 


1929 


1987 


1981 


1982 


1970 


3000 


2874 


2842 


2962 


2947 


2933 


2920 


4000 


3777 


3730 


3921 


3898 


3893 


3840 


5000 


4684 


4599 


4882 


4839 


4825 


4749 


6000 


5561 


5453 


5814 


5771 


5707 


5660 


7000 


6426 


6293 


6759 


6692 


6632 


6563 


8000 


7271 


7119 


7693 


7608 


7552 


7460 


9000 


8115 


7937 


8620 


8518 


8440 


8347 


10000 


8966 


8746 


9541 


9420 


9334 


9229 



"EVs r(-'>' gives the total number of different eigenvalues of T^^\ and "good 
EVs" means the number of different eigenvalues which are not spurious ac- 
cording to CuUum's and Willoughby's criterion. Note, however, that a "good" 
eigenvalue is generally only an approximation to an eigenvalue of but this 
approximation has not necessarily converged. 

A "good" eigenvalue turns into an "exact" one after convergence. Eig. 1 shows 
the convergence of the lowest positive eigenvalues and of the highest ones on 
an unquenched 164 lattice at /3 = 2.3, k, = 0.15. There exists an almost point 
reflection symmetry in that the picture for the highest negative and for the 
lowest eigenvalues looks practically the same if one changes signs at the axes. 
On other lattices the flgures look qualitatively similar. 

In the runs where complete spectra were determined (see Sec. 4) the following 
observation could be conflrmed [17]: Accurate approximations to eigenvalues 
of Q will be stabilized eigenvalues of T^^K This means that if an eigenvalue 
has converged for some j, then it will be an "exact" eigenvalue of T^-'^ for 
any larger value of j . Hence, we infer from Eig. 1 that not only the extremal 
eigenvalues converge fastest but also the eigenvalues close to zero. This fact 
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convergence of lowest positive A's 



# from 
lowest 
positive 




0.01 0.02 0.03 0.04 0.05 

A 



convergence of highest A's 



# from 
highest 




Fig. 1. Convergence of the lowest positive eigenvalues (top) and of the highest ones 
(bottom) of the unquenched Q as found by the Lanczos method on a IG"* lattice, 
(3 = 2.3, K = 0.15. From "o" to "•" the curves show the "good" eigenvalues of T^-') 
found with i = 500, 1000, 2 000, 3 000, 4 000, 5 000, 6 000, 7 000, 8 000, 9 000, and 
10 000 iterations, respectively. 
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was already pointed out in Refs. [13,14]. The reason is that the convergence is 
fastest in those regions of the spectrum where the spectral density is lowest, 
and in the configuration of Fig. 1 this is the case close to and furthest from 
the origin. We will make the connection between convergence properties and 
the local spectral density more precise in Sec. 4. 

In the regions of lowest spectral density the eigenvalues usually converge in an 
"ordered" fashion by which we mean the following. Between any two converged 
eigenvalues there is no "good" eigenvalue which has not yet converged also or 
which will not converge within a few more iterations. This observation is in 
accordance with Ref. [2] for staggered fermions. 

After an eigenvalue has converged, it starts to replicate [17] which means that 
it appears as a multiple eigenvalue of T^^\ even if it is only a simple eigenvalue 
of Q. The converse can also happen: A truly degenerate eigenvalue may only 
appear as a simple eigenvalue of T^-'^ for some j . In any case, the multiplicities 
of the eigenvalues of T^-'^ do not refiect the multiplicities of the eigenvalues 

of g. 

As mentioned in the introduction we are interested in a comparison of the 
partially converged Lanczos method with an accelerated conjugate gradient 
algorithm for low-lying eigenvalues [15,16]. This practical conjugate gradient 
variant determines low-lying eigenvalues of with rigorous practical error 
bounds. Analogous error bounds do not exist for the Lanczos method, ^ but 
we are able to monitor the accuracy of the non-spurious eigenvalues of T^-'^ 
a posteriori by referring to the converged (stabilized) values. We squared and 
sorted the "good" eigenvalues of T^-'^ and performed this convergence monitor- 
ing. By requiring a relative accuracy of 10""*, we obtain the results of Table 2. 
These results can be compared with Table 1 of Ref. [16], yielding an inferiority 
of the accelerated conjugate gradient method by a factor of 5 — 8 when one 
counts only the number of vector multiplications which are necessary for a 
relative error of 10""*. The iteration number j where this accuracy is reached 
is not known immediately, but we note that when the relative error is already 
10""*, then it will be reduced further very quickly. Hence, the entries in Table 2 
will not change much when we use the stabilization of the eigenvalues as an 
indication for their convergence. 

Instead of running the Lanczos procedure with A = Q one can alternatively 
use it with A = if one is not interested in the signs of the eigenvalues. 
Doing so, the results for the test configurations presented in Tables 1 and 2 
and in Fig. 1 remain similar. However, now one Lanczos iteration involves two 



^ CuUum and Willoughby quote error estimates [17], but they are not practical in 
lattice gauge theory, because either all eigenvalues of both T^^^ and r2"'^ are required 
or one needs eigenvector approximations. 
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Table 2 

Examples for the numbers of lowest eigenvalues of with a relative error < IQ-^^ 
found after j iterations of the Lanczos procedure for the operator Q. The lowest 
eigenvalue of is denoted by A^. 





6^ lattice 


63 • 12 lattice 


83 • 12 lattice 


84 lattice 




K = 0.15 


K = 0.15 


K = 0.15 


K = 0.20 




[3 = 1.80 


/3 = 2.12 


/3 = 2.12 


(3 = 0.00 




quenched 


unquenched 


unquenched 


quenched 


j 


\\ = 5.240 • 10-3 


Xl = 1.332 • 10-3 


Xl = 8.098-10-4 


Xl = 1.592 - 10-3 


250 














500 


7 


4 


1 





1000 


22 


21 


9 


2 


2000 


73 


56 


37 


14 


3000 


111 


85 


67 


28 


4000 


179 


145 


94 


38 


5000 


254 


182 


119 


50 






12^ lattice 


12^ lattice 


164 lattice 


164 lattice 




K = 0.15 


K = 0.15 


K = 0.15 


K = 0.15 




(3 = 1.80 


/3 = 2.40 


/3 = 2.12 


(3 = 2.30 




quenched 


quenched 


unquenched 


unquenched 


j 


= 4.752 • 10-3 


Xl = 1.364-10-4 


Xl = 7.703-10-4 


Xl = 1.367-10-4 


500 














1000 


3 


4 


2 


2 


2000 


16 


18 


8 


8 


3000 


25 


41 


17 


25 


4000 


37 


62 


23 


38 


5000 


49 


71 


30 


45 


6000 


57 


99 


49 


59 


7000 


70 


112 


63 


77 


8000 


79 


135 


75 


101 


9000 


101 


161 


83 


111 
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instead of one application of and the cost is roughly doubled. Therefore, it is 
advantageous to diagonalize the unsquared operator in case that one requires 
the lowest or highest eigenvalues of . This statement will become more and 
more significant, and the cost when working with will further grow, as a 
critical point is approached in the /3-/€-plane. This discussion is resumed in 
Sec. 6. 

As can be seen from Table 2, the number of converged eigenvalues is not a 
simple function of j . In particular there is no linear relationship between the 
number of converged eigenvalues and j, so that a proceeding as in Ref. [11] 
might require some refinement. A comparison of Tables 1 and 2 also yields no 
simple connection between the number of "good" eigenvalues and the number 
of converged eigenvalues. And finally, one can conclude from Table 2 that the 
computational cost for a certain fraction of the total number of eigenvalues 
cannot be expected to increase with less than the square of the lattice volume. 



4 Computation of complete spectra 

Although usually one does not want to compute all eigenvalues in production 
runs, it is nevertheless worthwhile to consider also complete spectra in the 
present technical study. One reason is that one can derive valuable information 
about the quality of the eigenvalues which are obtained for j <^ n. 

Let us start by looking at the convergence behavior of the Lanczos algorithm 
not only for the fastest converging eigenvalues but over the entire spectrum. 
We would like to concretize what is meant by the general statement of Ref. [17] 
that the local gap structure plays a role. An example with n = 20 736 eigen- 
values in an unquenched configuration on a 6"^ • 12 lattice at /3 = 2.12, k, = 0.15 
is shown in Fig. 2. The Lanczos method was run both with A = Q and with 
A = Q^. Every 500 iterations the matrices T^^^ and r2"'^ were diagonalized, 
and the "good" eigenvalues of T^-'^ were compared with reference data from 
2^(82944) ^ggg ^g^Q^ ^Yie choice of j = 82 944 = 32 x 6^ • 12). Every dot in 
Eig. 2 indicates a point where an eigenvalue A of Q or of has converged. 
Diamonds show the appropriately scaled spectral densities p (see Sec. 6). We 
notice that the number of Lanczos iterations required for convergence is in- 
timately related with the local spectral density. Eor instance, in the Lanczos 
procedure with A = Q we see a nice proportionality around A = 0. 

Apart from the example of Eig. 2, complete spectra were determined on a 
number of quenched and unquenched configurations with two fiavors of dy- 
namical fermions on 4"* - 8"^ • 12 lattices, and almost on a 12"* lattice. These 
configurations are quoted in Tables 4 and 5 below. In nontrivial gauge fields 
no degeneracy of any eigenvalue was found, neither for Q nor for Q^. In case of 
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convergence with Q versus spectral density of Q 
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convergence with versus spectral density of 

2.5 



1.5 



0.5 



<5- 



"ooocbooooooooc 



0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 

A2 



Fig. 2. Convergence of the Lanczos method versus spectral density for Q (top) 
and for (bottom) in an unquenched gauge field at /3 = 2.12, k = 0.15 on a 
6^ • 12 lattice: Dots give the number j of Lanczos iterations in multiples of n which 
are required for convergence of A or A^, and diamonds indicate the spectral densities 
scaled such that the maxima have the same value as the maxima of j'/ra. 
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the operator Q we always encountered an equal number of negative and posi- 
tive eigenvalues. Furthermore, no discrepancies were found when the squared 
eigenvalues of Q were compared with the results for . 

Examples for complete spectra are given in Fig. 3. The integrated densities 
of eigenvalues N{X) and N{\^) follow directly from the numerical data. They 
are normalized such that they take values between zero and one. Denoting on 
a lattice A of volume |A| the k-i\i (sorted) eigenvalue of Q by A/;, A; = 1, . . . , n, 
n = 4A^c|A|, then N is defined by 



and analogously for the squared operator. denotes the number of colours 
which is two in our case. 

We have the following consistency checks which provide good evidence that 
all computed complete spectra are correct. First, on all investigated 4"* - 8"^ • 12 
lattices the correct number of n different eigenvalues was found. Second, we 
have analytical sum rules for powers of the eigenvalues, cf. also Ref. [18]. The 
trace of Q equals zero, ^ and the trace of reads in any unitary gauge field 



Numerically we obtained |Tr Q|^10~®, and Tr came out with a relative ac- 
curacy of 10~® - 10~^^ (decreasing with increasing |A|). One could check for 
further sum rules by examining traces of higher powers of Q. However, except 
for Tr Q"^ = 0, compared with an absolute numerical value of <,10~'', such fur- 
ther checks were not performed. As another additional check for correctness 
one could compare the eigenvalues of T^-'^ either with those obtained by using 
a different initial Lanczos vector Vi or with those obtained in a gauge trans- 
formed configuration. Such checks were done for just one configuration and it 
turned out that the converged eigenvalues agree of course, but the spurious 
eigenvalues are generally dependent on Vi or they are gauge dependent. 

We confirmed that the numerical effort for the determination of complete 
spectra by means of the Lanczos algorithm grows with the square of the lattice 
volume. In order to obtain all eigenvalues in nontrivial gauge fields it was found 
empirically that j = 2n does not suffice, but j = An = 16A^c|A| worked in all 
cases, both for Q and for Q^. (Actually, the result of Fig. 2 suggests that 

^ Traces are understood over colour and spinor indices. 

In case of the squared staggered Dirac operator in SU(2) gauge fields, where 
every eigenvalue has a multiplicity of four, one obtains the complete spectrum when 
j equals only twice the number of different eigenvalues, i.e. when j = |A| [18]. 




(5) 





12 



integrated density Q 




integrated density 



iV(A2) 0.5 




Fig. 3. Integrated spectral densities of Q (top) and (bottom). The results are 
obtained in a quenched gauge field at /3 = 1.8 on a G"* lattice. At k = 0.15 curves of 
a quenched 12'^ gauge field at [i = 1.8 and unquenched gauge fields at [i = 2.12 on 
6"^ • 12 and 8"^ • 12 lattices coincide (on the scale of the figure) with the result on the 
G"* lattice. 
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j 2.5 n would also work.) With the choice of j = 4n the exact eigenvalues 
close to and furthest from the origin have multiplicities in T^-'^ of 0(10) - 
(9(100). However, the majority of the other correct eigenvalues of T^-'^ has 
a multiplicity of only one. Therefore it is generally too costly to accept an 
eigenvalue of T^-'^ as correct only when it has replicated. In the general case 
one better relies on the criterion that an eigenvalue has converged when it is 
not spurious and does not change with increasing j . 

Finally, we mentioned that the complete spectrum could almost be determined 
on a 12"* lattice. Because of memory limitations the largest matrix T^-'^ which 
could be handled is restricted to j = 21 • 12"* = 2.625 n. With this choice 
of j we investigated a quenched 12"* gauge field at /3 = 1.8, k, = 0.15. The 
computations required something like two weeks of workstation CPU time. 
(Six hours are required for the 100 lowest eigenvalues of Q^, cf. Table 1.) All 
eigenvalues of Q but four could be found, i.e. a total of 165884 out of 165888. 
For the traces these 165884 values yielded Xk = —0.09087 and A| = 
-0.04237 instead of zero, and E = 46612.04 instead of 46613.15. When 
one computes with these 165884 values the normalized integrated densities 
N[X) and A^(A^), one encounters no difference compared with the result on a 
6'^ lattice, only statistical noise which can be recognized on the 6'^ lattice is 
completely smoothed out on the 12"* lattice. 



5 Fermionic determinants 



In this section we give results for fermionic determinants of the configurations 
where all eigenvalues were computed. Besides the results for Wilson fermions 
we also quote values of determinants of staggered fermions whose spectra were 
determined in a previous publication [18]. It must be stressed that all numbers 
listed here are only examples for individual configurations. 



5.1 Wilson fermions 



We start by referring to free Wilson determinants in Table 3. In this case the 
eigenvalues of Q are known analytically. 



A = ± 



-j^ 4 4 

[7\ I]cosp^)2 + J^sin^p^ 

/i=l /i=l 



/(4 + -) (7, 



where each A appears twice with either sign, and = 2Trk^/L^ with 
0, 1, . . . , — 1 and is the extension of the lattice in //-direction. 
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Table 3 

Free Wilson fermions (in a trivial gauge field): (log;^Q det Q)/(4iVc|A|). In case ol 
K = 0.125 zero modes are excluded Irom the determinant. 



|A| 


K 


= 0.1 


K 


= 0.125 


K 


= 0.15 


K 


= 1/6 


K 


= 0.20 


K = 0.25 


44 


-0 


254927 


-0 


293846 


-0 


332655 


-0 


350106 


-0 


376619 


-0.639097 


6^ • 12 


-0 


253839 


-0 


295955 


-0 


329555 


-0 


346555 


-0 


371246 


-0.415416 


84 


-0 


253802 


-0 


295996 


-0 


329410 


-0 


346378 


-0 


370939 


-0.406257 


124 


-0 


253797 


-0 


296251 


-0 


329390 


-0 


346351 


-0 


370860 


-0.394140 


164 


-0 


253797 


-0 


296301 


-0 


329390 


-0 


346350 


-0 


370856 


-0.392127 


324 


-0 


253797 


-0 


296325 


-0 


329390 


-0 


346350 


-0 


370856 


-0.391263 



Table 4 

Examples for Wilson fermions in quenched SU(2) gauge fields: 

(log;^Q det Q)/(4iVc|A|). The case (3 = corresponds to a random gauge field. 



|A| 


/3 


K = 0.1 


K = 0.125 


K = 0.15 


K = 1/6 


K = 0.20 


K = 0.25 


64 


1.80 


-0.254757 


-0.299665 


-0.339289 


-0.362812 


-0.402336 


-0.440477 


64 


2.80 


-0.254331 


-0.298278 


-0.334999 


-0.354646 


-0.383828 


-0.409753 


64 


0.00 


-0.255279 


-0.301046 


-0.342455 


-0.368026 


-0.415065 


-0.477151 


84 


0.00 


-0.255277 


-0.301041 


-0.342445 


-0.368010 


-0.415038 


-0.476955 



Quenched results are contained in Table 4, while unquenched results are col- 
lected in Table 5. One concludes a nice exponential dependence of the deter- 
minant on the lattice volume, already for relatively small lattices. 

In a quenched Monte Carlo simulation the fermionic determinant is kept at 
a fixed value. From Table 5 one can get a feeling for the fiuctuations of the 
determinant in an unquenched run. We can compare the values of eight inde- 
pendent dynamical 63-12 configurations at /3 = 2.12, k = 0.15. The logarithmic 
entries in Table 5 fiuctuate by 0.0003 which translates to a fiuctuation of 
the determinant itself by six orders of magnitude on the 6^ • 12 lattice. 
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Table 5 

Examples for Wilson fermions in unquenched SU(2) gauge fields with two fiavours 
of dynamical fermions: (log;^Q det Q)/(4iVc|A|). 



|A| 


/3 


K 


= 0.15 


44 


1 


75 


-0 


338947 


6^ • 12 


2 


12 


-0 


337169 


6^ • 12 


2 


12 


-0 


337129 


6^ • 12 


2 


12 


-0 


337282 


6^ • 12 


2 


12 


-0 


337388 


6^ • 12 


2 


12 


-0 


337385 


6^ • 12 


2 


12 


-0 


337266 


6^ • 12 


2 


12 


-0 


337269 


6^ • 12 


2 


12 


-0 


337235 


8^ • 12 


2 


12 


-0 


337295 



Table 6 

Free staggered fermions (in a trivial gauge field): [log^o det(— ^^^^^^g + m^)]/(iVc|A|). 
In case of m = zero modes are excluded from the determinant. 



|A| 


171 = 


m = 0.05 


m = 0.2 


64 


0.878691 


0.835873 


0.853054 


124 


0.869726 


0.867211 


0.870592 


244 


0.868572 


0.868572 


0.871123 


324 


0.868501 


0.868615 


0.871133 


484 


0.868471 


0.868628 


0.871135 


644 


0.868465 


0.868630 


0.871135 


804 


0.868463 


0.868630 


0.871135 
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Table 7 

Examples for staggered fermions in nontrivial SU(2) gauge fields: 

[log;^Q det(— _^gj^g + m^)]/(iVc|A|). The case (3 = corresponds to a random gauge 

field. 



|A| 


/3 


171 = 


m = 0.05 


m = 0.2 


quenched 


D 


0.0 





498590 


0.511586 





553658 


D 


"1 o 

1.8 





639018 


0.645632 





670574 


1^4: 

D 


2.8 





785223 


0.785556 





-7 ri A o o '1 

790327 


12 


0.0 





499169 


0.513180 





r r A cxcw 

554995 


124 


1.8 





640427 


0.648199 





672902 


124 


2.0 





674723 


0.680302 





699873 


124 


2.4 





757622 


0.758471 





765383 


124 


2.6 





777512 


0.777890 





782955 


124 


2.7 





784554 


0.784866 





789400 


124 


2.8 





788623 


0.788923 





793297 


unquenched 


64 


1.8 






0.697371 


0.704457 


64 


2.8 






0.799398 


0.802219 


124 


1.8 






0.638295 


0.707586 


124 


2.4 






0.717424 


0.776017 


124 


2.8 






0.779199 


0.798827 



5.2 Staggered fermions 

The Lanczos algorithm was studied earlier for the staggered Dirac operator 
Pstag (also with periodic boundary conditions) in Ref. [18], but no determi- 
nants were computed there. Here we use the old data and for a sake of com- 
pleteness quote results in Tables 6 and 7. Note that the staggered quark mass 
m is measured in units of twice the separation of neighbouring sites,® i.e. in 
a more frequent convention one would quote m = 0, 0.025, 0.01 instead of 0, 
0.05, 0.02, respectively. 

® The reason is that free staggered fermions enjoy translation invariance only by 
shifts of two lattice spacings [20]. 
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We add as a remark that the staggered data of Ref. [18] were further analysed 
by Halasz and Verbaarschot in [8]. These authors showed that the eigenvalue 
correlations of the staggered Dirac operator in SU(2) gauge fields are very well 
described by the Gaussian symplectic ensemble of a random matrix model with 
the chiral symmetry of the Dirac operator built in. 



6 Spectral densities 

The density of eigenvalues of around zero can be related with the chi- 
ral limit, if one connects the spectral density of with that of the Dirac 
operator (not multiplied by 75) in the spirit of the "pion norm" of Ref. [1]. 
Such an analysis will however be done elsewhere. From the algorithmic point 
of view spectral densities are interesting because they determine convergence 
properties, as we saw from Fig. 2. 

Figs. 4-6 show results for normalized spectral densities as obtained in in- 
dividual gauge field configurations of the present algorithmic investigation. 
The data points in these figures were obtained by a simple binning procedure, 
i.e. approximating p{\) = dN(A)/dA by AiV(A)/AA with AA = 0.01, and 
analogously for the spectral density p{\^) of . The solid curves indicate ev- 
erywhere reference data of unquenched gauge fields at n = 0.15, fi = 2.12 on 
6"^ • 12 and 8"^ • 12 lattices, where numerical results coincide. 

The information in each of Figs. 4 - 6 is redundant because of the relation 
p{\^) = /j(A)/(2A), which can be verified numerically to a very high precision. 
Despite this redundancy we consider it worthwhile to present results for p{X) 
as well as for /o(A^). Note for instance that in the K-range investigated p{\^) in 
Fig. 4 stays finite as A^ — )• because p{\ = 0) = 0. On the other hand, p{\^) 
must diverge for A^ — )• if p{\ = 0) is finite. This is the case for n = 0.25 in 
Figs. 5 and 6, and the case n = 0.20 in Fig. 5 seems to be close to a marginal 
situation. 

Remarkably, wherever we are able to compare data from different lattice sizes 
> G"* which correspond to the same physical situation, we cannot find finite 
size effects in the normalized spectral densities (using AA = 0.01), and also in 
the integrated densities on the scale of Fig. 3; the only effect of a larger lattice 
is to smooth out statistical fiuctuations. This is true for the unquenched gauge 
fields at n = 0.15, fi = 2.12 on 6"^ • 12 and 8"^ • 12 lattices, for quenched gauge 
fields at fi = 1.8, /€ = 0.15on6'* and 12"* lattices, and for random gauge fields 
(/3 = 0) on G"* and 8"* lattices in the K-range 0.10, . . . , 0.25. 
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density Q, /3 = 2.8 




density , /3 = 2.8 




Fig. 4. Spectral density of Q (top) and of (bottom) as obtained in a quenched 
gauge field at /3 = 2.8 on a 6"^ lattice. The solid lines indicate the densities in 
unquenched gauge fields at k = 0.15, (3 = 2.12 on 6"^ • 12 and 8"^ • 12 lattices. 
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density Q, /3 = 1.8 




-1 -0.8 -0.6 -0.4 -0.2 0.2 0.4 0.6 0.8 1 

A 



density Q'^ , /3 = 1.8 




Fig. 5. Spectral density of Q (top) and of (bottom) as obtained in a quenched 
gauge field at /3 = 1.8 on a G"* lattice. The curves at k = 0.15 coincide with the 
results on a quenched 12'^ lattice at /3 = 1.8 (except that fluctuations are smoothed 
out). The solid lines indicate the densities in unquenched gauge fields at k = 0.15, 
13 = 2.12 on 6^ • 12 and 8^ • 12 lattices. 
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density Q, /3 = 0.0 




density , /3 = 0.0 





K 


= 0.25 


+ 


H 


= 0.20 


X 


H 


= 1/6 


O 


H 


= 0.15 


□ 


H 


= 0.125 


A 


H 


= 0.10 




Fig. 6. Spectral density of Q (top) and of (bottom; note the different scale 
of the ordinate compared with Figs. 4 and 5) as obtained in a random gauge field 
(P = 0.0). There is no distinction between 6"^ and 8"* lattices, except that fiuctuations 
are smoothed out on the larger lattice. The solid lines indicate the densities in 
unquenched gauge fields at k = 0.15, /3 = 2.12 on 6"^ • 12 and 8"^ • 12 lattices. 
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Finally, let us resume the discussion about the relation between convergence 
and the spectral density. We noted that the complete spectrum could be ob- 
tained for any of the configurations studied if j = 4n, both for Q and for . 
This means that a complete computation is not affected by the degree of crit- 
icality of the system, because we included random gauge fields with n = 0.25 
which is believed to be a critical point (in strong coupling and in the large 
limit) [2f]. We also studied random configurations with n in the neighbour- 
hood of 0.25 and up to 2.0, and we did not encounter any exceptional example 
where not all eigenvalues were detected with j = An. However, the picture 
changes if we pay attention only to some of the low-lying eigenvalues of Q^. 
In this case one can witness quite a drastic difference between the Lanczos 
procedure with A = Q and A = in the vicinity of a critical point. It takes 
much longer until the first few eigenvalues have converged when one runs the 
Lanczos method with instead of Q. The reason is obvious when one looks 
at the spectral densities near zero in Fig. 6 and recalls Fig. 2. 



7 Conclusions and comparison with the algorithm of Ref. [16] 

If one is willing to pay the necessary CPU time — which increases from 2.5 
workstation minutes on a 4"* lattice to more than two workstation weeks on 
a 12"* lattice — one can compute the complete spectrum of Q in quenched 
and unquenched SU(2) gauge fields by means of CuUum's and Willoughby's 
Lanczos procedure. Such a computation of the complete spectrum is also reli- 
able in cases where the configuration is in the immediate neighbourhood of a 
critical point, as was seen in random gauge fields around k, = 0.25. We always 
detected the correct number of eigenvalues which equals the dimension of the 
Dirac operator (except for four missing eigenvalues on the 12"* lattice, but this 
was due to memory limitations), and we were able to check the correctness of 
these eigenvalues by analytical sum rules. 

With data for complete spectra one can compute fermionic determinants and 
spectral densities. For the determinants we found a nice exponential depen- 
dence on the lattice volume already for relatively small lattices. The spectral 
densities are expected to contain information about the question of chiral sym- 
metry breaking and the phase structure of the theory. Possibly, with data for 
all eigenvalues one can also investigate questions of universality similarly to 
Ref. [8]. However, such physical applications will be addressed elsewhere. 

In Fig. 2 we showed that the convergence behaviour of the Lanczos method is 
intimately linked with the local spectral density p over the entire spectrum. 
Convergence is fastest in those regions of the spectrum where the density of 
eigenvalues is lowest. Thus, if one is interested in the low-lying eigenvalues of 
one should run the Lanczos procedure with Q and not with Q^, because 
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of the relation p{X^) = p{X) / {2X). The difference in performance is roughly a 
factor of two in the computational cost (measured in the number of vector 
multiplications) when the system studied is not very close to criticality, oth- 
erwise the cost ratio gets even bigger. 

Finally, let us turn to a comparison of the Lanczos procedure with the ac- 
celerated conjugate gradient algorithm of Ref. [f6] for the computation of 
low-lying eigenvalues of Q^. In Ref. [16] the same unquenched configurations 
as in the present paper were studied and a number of low-lying eigenvalues 
with a relative accuracy of 10"'* were computed. We find the following results. 
First, by means of the reference data of [16] we are able to confirm that in 
the tests performed the Lanczos procedure with a small number of iterations 
j <^ n never skipped an eigenvalue with small modulus. The reason is that in 
the regions of lowest spectral density the Lanczos method shows an "ordered" 
convergence by which we mean that between any two converged eigenvalues 
there is no eigenvalue which has not yet converged also or which will not 
converge within a few more iterations. Second, comparing the performance 
measured in vector multiplications (which require the major fraction of 
the total CPU time), the accelerated conjugate gradient algorithm is inferior 
to the Lanczos method by about a factor of 5 — 8. Possibly, this factor might be 
larger for systems close to criticality where the spectral density favours com- 
putations with Q instead of Q^. Hence, a practitioner may prefer the Lanczos 
method provided he has a priori information about degeneracies in the spec- 
trum. However, from a rigorous point of view the conjugate gradient approach 
is favourable because it yields all eigenvalues with a rigorous straightforward 
and practical error bound, and all eigenvalues are detected with their correct 
multiplicities. Furthermore, approximations to eigenvectors are obtained as a 
by-product. This latter point is a clear advantage for applications where not 
only eigenvalues but also eigenvectors are needed, e.g. when one is interested 
in the contribution of the low-lying eigenmodes to physical observables. 
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